import os

def get_md_density():
    os.system("readener.py Density > density.log")
    molecular_group = {}
    molecular_md_density = {}
    f = open("density.log", 'r')
    for i in f:
        a = i.split()
        if len(a) > 1:
            b = a[0].split("/")
            if len(b) > 2:
                molecular_group[b[2]] = b[1]
                molecular_md_density[b[2]] = float(a[3]) / 1000
    return molecular_group, molecular_md_density

def get_exp_density():
    molecular_exp_density = {}
    f = open("exp_data.log", 'r')
    for i in f:
        if i.strip().startswith("#"):
            pass
        else:
            a = i.split()
            print a
            molecular_exp_density[a[1]] = float(a[4])
    return molecular_exp_density

def write_density_file():
    molecular_group, molecular_md_density = get_md_density()
    molecular_exp_density = get_exp_density()

    o = open("density_analysis.log", 'w')
    o.write("%20s%20s%8s%8s%8s%8s\n"%("#group", "molecule", "MD", "EXP", "dev", "abs_dev"))
    for i in molecular_group.keys():
        o.write("%20s%20s"%(molecular_group[i], i))
        o.write("%8.3f%8.3f"%(molecular_md_density[i], molecular_exp_density[i]))
        dev = (molecular_md_density[i] - molecular_exp_density[i]) / molecular_exp_density[i]
        abs_dev = abs(dev)
        o.write("%8.3f%8.3f\n"%(dev, abs_dev))
    o.close()

if __name__ == "__main__":
    write_density_file()
